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Excitons, trions, biexcitons, and exciton-trion complexes in two-dimensional transition metal dichalcogenide 
sheets of M 0 S 2 , MoSe 2 , MoTe 2 , WS 2 and WSe 2 are studied by means of density functional theory and path 
integral Monte Carlo method in order to accurately account for the particle-particle correlations. In addition, 
the effect of dielectric environment on the properties of these exciton complexes is studied by modifying the 
effective interaction potential between particles. Calculated exciton and trion binding energies are consistent 
with previous experimental and computational studies, and larger systems such as biexciton and exciton-trion 
complex are found highly stable. Binding energies of biexcitons are similar or higher than those of trions, 
but the binding energy of the trion depends significantly stronger on the dielectric environment than that of 
biexciton. Therefore, as a function of an increasing dielectric constant of the environment the exciton-trion 
complex ’’dissociates” to a biexciton rather than to an exciton and a trion. 


Layered transition metal dichalcogenides (TMD) are chem¬ 
ically, thermally, and mechanically stable even in the mono- 
layer form, and thus, provide an ideal platform for study¬ 
ing condensed matter physics in two dimensions. The semi¬ 
conducting TMDs present many unusual optical properties 
such as strong excitonic effect^, valley-dependent circular 
dichroisrrP, and second-harmonic generatiorP, whose mag¬ 
nitude depends sensitively on the number of layers. For in¬ 
stance, the prototypical M 0 S 2 material is a semiconductor 
with 1.1 eV indirect band gap in bulk, but 1.9 eV direct band 
gap in the monolayeiEl. Importantly, the reduced dimensional¬ 
ity is manifested in a large exciton binding energy of 0.5-1 eV, 
but also of significant binding energy in the case of charged 
excitons, or trions, consisting of three charge carriers. This 
suggests that even larger complexes might be stable. Indeed, 
first experimental reports assigned to biexciton formation have 
very recently appeared in the literaturd^. 

Theoretical studies are invaluable in predicting the stabil¬ 
ity of these complexes and in interpreting the experimen¬ 
tal results. Excitons can be calculated reliably from first- 
principles by solving the Bethe-Salpeter equation (BSE) on 
top of quasi-particle band structure. Binding energies have 
also been calculated using simple variational or tight-binding 
models based on an effective 2D interaction potential and the 
effective mass approaclPtlll yielding fairly good agreement 
with experiments and with the BSE results in the case of exci¬ 
tons. This has also raised interest t o ap ply similar approaches 
to study larger exciton complexeg^ES in comparis on to the 
theoretical estimates based on quantum well systemJ^EIIll 

Difficulties in constructing reasonable wave function ansatz 
in the case of the larger complexes hinders straightforward 
extension of the simple variational models. Within the effec¬ 
tive mass approach, quantum Monte Carlo (QMC) methods, 
such as diffusion Monte Carlo and path integral Monte Carlo 
(PIMC), provide accurate and powerful means for studying 
few-particle systemJ^^ES. The main advantage in QMC meth¬ 
ods is the exact account of particle-particle correlations, which 
is particularly important in accurate description of exciton 
complexes. To this end, we utilize the PIMC method, which is 
a basis set free approach for solving finite temperature quan¬ 


tum statistics. 

Although 2D materials are not directly bonded with the 
environment, due to their thinness, they are highly sensitive 
to electromagnetic fields, doping, or dielectric screening of 
their surroundings. In particular, the Coulomb interaction 
between an electron and a hole in exciton is screened by 
the dielectric envir onment and the binding energy changes 
dramatic alljUSHiZlz^Hl p j-ather surprising then that the 
effect of dielectric environment on the binding energy of 
trions, let alone on the larger complexes, has been rarely 
investigateci^. 

Here, we present the results from PIMC simulations for ex¬ 
citon, trion, biexciton, and exciton-trion complexes for a set 
of the most common layered TMD materials: M 0 S 2 , MoSe 2 , 
MoTe 2 , WS 2 , and WSe 2 . We focus on the binding energies 
and mean particle distances. Our approach is based on the 
effective 2D interaction and effective masses, for which rele¬ 
vant parameters are calculated using density functional theory 
(DPT). In addition, the effect of environment is accounted for 
in the interaction potential, which allows us to demonstrate the 
effect of the surroundings to the binding energies and mean 
distances. 

The geometry of the system is schematically shown in Pig. 
[2 The electrons and holes are confined to the 2D sheet placed 
on top of a substrate. We consider two-dimensional many- 
body Hamiltonians where the interaction potential between 
charged particles is given as in Refs. l9llOI23l 



where Hq and Fq are the Struve function and the Bessel func¬ 
tion of the second kind, respectively, qi represents the charge 
of the particle (here qi = e or qi = —e). The length scale in 
the presence of dielectric environment is given as 

r* = 27rx2D/«:- (2) 

X2D is the 2D polarizability of the sheeP®, which can be eval¬ 
uated from the bulk in-plane dielectric constant £|| and the 
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FIG. 1; (a) Schematic illustration of the system geometry, (b) Illus¬ 
tration of the exciton, trion, and biexciton systems constructed from 
the electrons on the conduction band and holes in the valence band, 
(c) Effective interaction potential in the case of Senv = 1, 3.9 (silica), 
and 10 (sapphire), as compared to Coulomb interaction in homoge¬ 
neous system scaled by k. 


layer separation ds (here half of the perpendicular lattice con¬ 
stant since all considered materials have two layers within the 
primitive cell); X 2 D = '^s(£|| — l)/47r. k is the average di¬ 
electric constant of the environment. Here we consider the 
experimentally most relevant case of TMD sheet placed on 
a substrate of dielectric constant Eenv and vacuum or air on 
the other side, in which case k = (1 + eenv)/2. The inter¬ 
action potentials for three representative values of k in Fig. 
[^c) show strong screening by the 2D sheet at short distances 
and approaching 1/kt at r > r*. We note, that this type of 
potential is only valid when k < 

The four considered exciton complexes are illustrated in 
Fig.[TJb). In monolayer TMDs, due to the spin-orbit coupling 
and the lack of inversion symmetry, the valence band maxi¬ 
mum (VBM) at K-point is split to two bands with distinct spin 
orientation. The spins are further coupled to the valley index 
(K or K’). Thus, if the many-body system is supposed to have 
two holes of opposite spins in the topmost valence band state, 
they must be located in different valleys. The conduction band 
minimum (CBM) is nearly degenerate (here assumed degen¬ 
erate) and there is no such restrictions for the electrons. 

For the model interaction potential two material specific pa¬ 
rameters are then needed; polarizability of the sheet and the 
effective masses of electrons and holes. Both of these quanti¬ 
ties can be obtained quite reliably from DFT calculations. Our 
calculated numbers are collected in Table |I] The atomic struc¬ 
tures are optimized using revB86b-DF2 functional^, which 


yields structural parameters in very good agreement with the 
experimenP^. Using these structures, the effective masses and 
dielectric constants are then calculated using the PBE func¬ 
tional. Bulk dielectric constants are calculated using density 
functional perturbation theory. Since all materials considered 
here have direct gap with VBM and CBM located at K-point 
(in the monolayer form), the effective masses are obtained 
from the monolayer structure by fitting parabolas to the K- 
valley and accounting for the spin-orbit couplin^^. 


TABLE I; Material parameters needed for the effective mass theory 
description of the exciton complexes. Layer separation ds, and X 2 d 
are in A. The two values for rrih correspond to the two spin-orbit 
split valence bands at K (first for the highest band), e is dielectric 
constant in directions parallel to the plane and normal to the plane. 


Material 

ds 

rrih 

rris 

£ 

X2D 

M 0 S 2 

6.180 

-0.54/-0.61 

0.47 

15.46/6.46 

7.112 

MoSe2 

6.527 

-0.59/-0.69 

0.55 

17.29/7.95 

8.461 

MoTe2 

7.054 

-0.62/-0.75 

0.57 

21.87/11.02 

11.715 

WS 2 

6.219 

-0.35/-0.49 

0.32 

13.92/5.98 

6.393 

WSe2 

6.575 

-0.36/-0.53 

0.34 

15.47/7.22 

7.571 


The path integral Monte Carlo simulations are carried out at 
T = 10 K using the effective mass approach with the masses 
obtained from our DFT calculations, and the effective inter¬ 
action potential given in Eq. Q- Apart from the exciton- 
trion complex we can consider our particles as "boltzman- 
nons", i.e., they obey the Boltzmann statistics and are treated 
as distinguishable particles. In case of negatively charged tri- 
ons this is possible by assigning spin-up to one electron and 
spin-down to the other one. With biexcitons we simply ap¬ 
ply the same for the positive particles, also. The exciton-trion 
complex requires the account of Eermi statistics, which in this 
work is dealt with by the restricted path integral Monte Carlo 
approach, and the free particle nodal restrictiorP^. 

In this work the statistical standard error of the mean with 
2(7 limits is used as an error estimate for all observables from 
our PIMC simulations. Sampling in the configuration space 
is carried out using the Metropolis procedure's with multi¬ 
level bisection move^, and the thermal estimato^Hl is used 
in the calculation of the total energy. We employ the PIMC 
method with the primitive approximatioiP^, for which we find 
that T = 10 K describes the ground states of our systems 
accurately, and that using Trotter number M = 4000 yields 
good balance between feasible amount of computer time and 
accuracy. As the Trotter number M tends to infinity, the ex¬ 
act many-body results are obtained, but high-accuracy is often 
found with reasonable values of M. More details are given in 
the Supplemental Material (SMj^. 

Before discussing the binding energies, it is useful to illus¬ 
trate the spatial distribution of the electrons and holes in these 
systems. In Pig. we show the ^-coordinate distribution of 
all particles, when x-axis is chosen to pass through two of the 
particles with origin at their center-of-mass. In exciton, the 
electron distribution when hole is fixed at the origin is sim¬ 
ilar to that of hydrogen Is state, showing exponential decay 
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at larger distance. In the case of negative trion, the hole is 
largely located between the electrons. When the distance be¬ 
tween electrons in biexciton is large, there is one hole located 
close to each electron (solid lines in Fig. |^). When the dis¬ 
tance between electrons is small (dashed lines in Fig.[^), the 
extent of the hole wave functions becomes too large to make 
such distinction any more. 

The root-mean square (rms) averaged electron-hole separa¬ 
tion for exciton in suspended M 0 S 2 is found to be 11 A, which 
is in agreement with previous model calculation and with 
the GVFh-BSE calculated rms radius of 1 nm reported in Ref. 
[3n The electron-electron rms separation for A“ is 29 A and 
for biexciton 23 A (see SM for tabulated data of inter-particle 
distances in all the considered tmdP). 



x(A) 



-40 -20 0 20 40 

x(A) 


FIG. 2: Particle coordinate distributions in the case of (a) exciton, 
(b) trion, and (c) biexciton. In (a), the hole is fixed at the origin. 
In (b) and (c), x-axis is chosen to pass along two particles as illus¬ 
trated. Origin is at the center-of-mass of these two particles. In (c), 
the distributions for electron-electron distance larger (solid lines) and 
smaller (dashed lines) than 40 bohr are distinguished. 


The binding energies for all complexes are given in Table [H] 
together with the experimental results. Comparing the num¬ 
bers over different materials, we observe that the binding en¬ 
ergy depends strongly on the material polarizability X 2 D, but 
weakly on the effective masses. As a consequence, the results 
for A and B cases are generally very similar. For instance, 
even if the effective mass ratio in WSe 2 differs by more than 
30%, the binding energy differs less than 5%. Our calculated 
results for B~ and 5+ trions and BB biexcitons were also 
very similar to the A counterpart and thus omitted from Ta¬ 
ble |I^ Finally, the effective mass insensitivity also extends to 


negative and positive trions having nearly identical binding 
energies. 

When comparing to the results obtained using ah initio 
GVFh-BSE approach, as collected to Table [n| our model tends 
to underestimate the exciton binding energies by 0.1-0.2 eV; 
a satisfactory agreement. The difference can arise from the 
fairly small extent of the exciton wave function, and con¬ 
sequently breakdown of the effective mass approximation. 
There are no first principles results available for comparison 
in the case of trion, biexciton, or exciton-trion complex. Nev¬ 
ertheless, our results for excitons and trions are close to those 
obtained using the variational model for the case of k = 
but clearly smaller than the estimate for biexciton given in 
ReflU 

In experiments, TMD sheets are rarely suspended in vac¬ 
uum. Dielectric environment has strong effect on the screen¬ 
ing of the interactions within the sheet and consequently on 
the binding energieJ ^^^ I ^' l ^^ l The binding energies as a func¬ 
tion of the average dielectric constant of the environment k 
are shown in Eig.|^ Eor excitons, the dependence on k can be 
fitted reasonably well with Ei,{k = when a « 0.7 in 

agreement with the asymptotic form found in Ref. With 
increasing k, the trion binding energy decreases faster than 
that of biexciton and they are found equal at k « 4. Since 
the short range interaction is more strongly affected by the di¬ 
electric environment (cf. Eig. [^), and considering the trion 
geometry shown in Eig.|^, the repulsive electron-electron in¬ 
teraction is affected less than the attractive electron-hole in¬ 
teraction. If biexciton is approximated as two weakly bound 
excitons, then their binding energy should have only weak k- 
dependence. 




FIG. 3: (a) Calculated exciton binding energies as a function of 

K together with error bars (2a). Dotted line is guide for the eye. 
Two fits are also shown: The coefficients for the first (solid line) are 
a = —325 meV and b — 853 meV. The second fit (dashed line) was 
proposed in Ref. (b) Binding energies and error bars for trion, 
biexciton, and exciton-trion complex as a function of k. 

Comparing our calculations to experimental results is hin¬ 
dered by the large variations in the reported numbers. For¬ 
tunately, however, many of these experiments are carried out 
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TABLE II: Binding energies for all considered systems both in vacuum and in dielectric environment described by k = 2. Experimental and 
computational results from literature are collected for comparison. Exciton energies are given in eV and binding energies of exciton complexes 
are given in meV. The 2a statistical error estimate is given in parentheses for the PIMC results. 



M 0 S 2 

MoSe2 

MoTe2 

WS 2 

WSea 

Exciton A 

0.5265(2) 

0.4769(2) 

0.3752(2) 

0.5098(2) 

0.4564(2) 

Exciton B 

0.5339(2) 

0.4853(2) 

0.3828(2) 

0.5309(2) 

0.4777(2) 

Other calc. 

o.spl, o.'P, o.spl 

O.bP 


O.P 


Exciton A at K = 2 

0.3486(2) 

0.3229(2) 

0.2608(2) 

0.3229(2) 

0.2946(2) 

Expt./Other calc. 

O.4PIII, 0 . 4 pl 


O.5P 

0.3P 

O. 3 P 

Trion A“ 

32.0(3) 

27.7(3) 

21.0(2) 

33.1(3) 

28.5(3) 

Trion A+ 

31.6(3) 

27.8(3) 

20.9(3) 

33.5(4) 

28.5(4) 

Trion A~ at k = 2 

24.7(3) 

22.1(3) 

17.1(2) 

24.3(3) 

21.5(3) 

Expt. 

iP 

2^401421 

2P 

3P 

3]E1E3 

Biexciton A A 

22.7(5) 

19.3(5) 

14.4(4) 

23.9(5) 

20.7(5) 

Biexciton A A at k = 2 

20.3(5) 

17.4(4) 

12.9(4) 

20.9(5) 

18.1(4) 

Ex-l-Trion AA~ 

17.0(6) 

16.4(5) 

12.5(5) 

14.9(6) 

14.9(6) 

Ex-l-Trion AA“ at k = 2 

13.5(4) 

12.7(4) 

10.0(4) 

13.1(5) 

12.2(4) 


with the TMD sheet placed on a Si 02 substrate. Then, in order 
to facilitate the comparison we have recalculated all binding 
energies for k = 2, approximately corresponding to a situa¬ 
tion of Si 02 on one side and vacuum/air on the other. The 
results are given in Table [n| In the case of excitons, our cal¬ 
culated binding energies are now within 0.1 eV of the exper¬ 
imental values. For trions, our calculations are also generally 
in line with the experimental results or slightly underestimated 
as in the case of excitons. It is worth noting though, that 
some authors have reported similar binding energies for the 
A” and A+ trionJSl in agreement with our results, whereas 
others have found larger binding energy for the A^ than for 

the A+EiE5l. 


have similar electron-hole exchange contributions as excitons. 
However, if we consider that the electron at K is bound to hole 
at K’ and vice versa, these are dark excitons with no electron- 
hole exchange. Comparing such conhguration to the energies 
of two bright excitons that are missing exchange energies in 
our calculations, leads to total energy correction of 40 meV. 
Thus obtained binding energy of 60 meV is then in line with 
the experimental values. Another explanation could be that in 
the experiment the biexcitons are bound to e.g. impurities. 

The only report of exciton-trion complex, to the best of our 
knowledge, is in Ref. |4T] They deduced binding energy of 
4 ± 1.5 meV for MoSe 2 on Si 02 . Our number from Table [II| 
for this case is 12.7(4) meV, which is substantially larger. 


In contrast to excitons and trions, which are commonly 
observed, the reports for biexcitons are very scarce. You et 
al. report binding energy of 54 meV for both inter- and in¬ 
travalley biexcitons in WSeJ^, Mai et al. report binding en¬ 
ergy of 70 meV for AA’ biexciton (where A’ denotes exciton 
in the K’ valley) in M 0 S 5 R and Sie et al. found binding en¬ 
ergy of 40 m eV for AB biexciton, 60 meV for AA’ biexciton 
in MoSjP^. These experimental results are clearly larger than 
our calculated values of about 20 meV. 

Likely explanation to the discrepancy is the neglect of ex¬ 
change, and especially the electron-hole exchange, in our cal¬ 
culations. GIFh-BS E cal culations yield dark exciton 20 meV 
below the bright ond^^HI which originates from electron-hole 
exchange due to vanishing splitting of the conduction band in 
MoS 2 . This is the exchange energy that should be added to our 
calculated value to yield the bright exciton energy. Neglect of 
exchange between electrons in negative trion, i.e., setting their 
spins the opposite, inevitably leads to missing the exchange 
from one electron-hole pair, with energy contribution similar 
to an exciton or smaller since the electron and hole are more 
separated. In AA’ biexciton there is no exchange between 
electrons or between holes, since they are assumed to have op¬ 
posite spins. The electron-hole pairs at K or K’ valley should 






environment dielectric constant k 


FIG. 4: Selected average inter-particle distances for exciton (a), 
trion (b), biexciton (c), and exciton-trion complex (d) in M 0 S 2 as a 
function of k. 

Finally, to illustrate the dependence of the system size on 
the environment, the inter-particle distances as a function of 
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K are shown in Fig. The particle distances generally scale 
linearly with k. Trion and exciton-trion complex are found to 
dissociate after k > 8. We note, that the dissociation is facil¬ 
itated by the non-zero temperature (10 K) used in our PIMC 
calculations. Interestingly, the exciton-trion complex does not 
dissociate to exciton and trion, but to biexciton (compare to 
hi — /i 2 distance in biexciton) and a free electron. This is 
further confirmed by inspection of other particle-particle dis¬ 
tances, i.e., electron-hole distances (not shown in the figure). 

In conclusion, we have studied excitons, trions, biexci- 
tons, and exciton-trion complexes in two-dimensional transi¬ 
tion metal dichalcogenides by an approach combining density 
functional theory with quantum Monte Carlo method. We fo¬ 
cused on the binding energies, inter-particle separations, and 
on the role of dielectric environment. Our approach repro¬ 
duced exciton and trion properties in reasonably good agree¬ 


ment with experiment. We found that the larger complexes 
should also be stable with binding energies comparable to 
those of trions, although the relative stability can be controlled 
by the dielectric environment of the 2D sheet. Due to the 
large binding energies, environmental control, and coupling 
with the valley and spin indices of the material, we expect 
TMD materials to provide a versatile “laboratory” for study¬ 
ing, experimentally and theoretically, the physics of correlated 
many-body systems going even beyond the 3-5 particle com¬ 
plexes considered here. 
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through its Centres of Excellence Programme (2012-2017) 
under Project No. 251748 (HPK). We also thank CSC-IT 
Center for Science Ltd. and Tampere Center for Scientific 
Computing for generous grants of computer time. 
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